Julia数据科学系列-StatsBase包

datascience, doc, julia, statistics    
  1. 加权向量AbstractWeights
    1. 加权向量AbstractWeights的子类型
    2. 支持方法
    3. 构造函数
  2. 标量统计
    1. 加权求和求平均
    2. 进阶平均值计算
    3. 矩/动差(moments)和积累量(cumulants)
    4. 积累量:
    5. 变异
    6. Z-score
    7. 熵相关
    8. 分位数相关
    9. 众数
    10. 概括统计量
    11. 信度分析
    12. 稳健性统计
    13. 偏差统计
  3. 散布矩阵和协方差

原文戳我

加权向量AbstractWeights

StatsBase引入AbstractWeights类型用来表述加权向量, 有两个好处:

  1. 与其他向量区分开, 可以更方便地配置函数传参

  2. 涉及权重的统计, 经常需要进行加权, AbstractWeights类型中包含加权信息, 不用重复计算

Note

  • 加权后的向量是对原始向量的封装, 并不拷贝原始向量

  • 加权信息是在构造时就计算好的, 所以如果加权信息已知, 可以在构造的时候就当作参数进行传递: AnalyticWeights(vs, wsum=sum(vs))

加权向量AbstractWeights的子类型

  • AnalyticWeights: 通常是[0-1]之间的值, 表示每个观测值的重要性, 构造函数简写aweights

  • FrequencyWeights: 表示每个观测值出现的次数, 构造函数简写fweights

  • ProbabilityWeights: 也叫采样权重, 表示对采样方法的矫正, 构造函数简写pweights

  • UnitWeights: 所有权重都是1, 等同于没有权重的向量, 用于保持输出类型稳定吧我猜, 构造函数简写uweights

  • Weights: 通用的权重类型, 构造函数简写weights

  • Exponential weights: 指数型增减的权重 eweights

支持方法

eltype
length
isempty
values
sum

构造函数

# AnalyticsWeights, FrequencyWeights, ProbabilityWeights, Weights都一样:
XXXWeights(vs, wsum=sum(vs))

UnitWeights{T}(s) # 构造T类型的长度为s的单元权重

# 构造函数的配套简写方法, 不支持显式传入wsum:
aweights(vs) # fweights, pweights, weights一样
uweights(s::Integer) # s => 长度, 构建Int类型的uweights
uweights(::Type{T}, s::Integer) where T<:Real

eweights复杂一些:

eweights(t::AbstractVector{<:Integer}, λ::Real; scale=false)
eweights(t::AbstractVector{T}, r::StepRange{T}, λ::Real; scale=false) where T
eweights(n::Integer, λ::Real; scale=false)

具体n, λ的意义暂略, 用到的时候再学吧!

标量统计

加权求和求平均

扩展基础统计函数以支持AbstractWeights:

# Base.sum, Base.sum!
sum(v::AbstractArray, w::AbstractWeights{<:Real}; [dims])
sum!(R::AbstractArray, A::AbstractArray, w::AbstractWeights{<:Real}, dim::Int; init::Bool=true)

# 新增 wsum 和 wsum! 函数
# ?? 既然都扩展sum了 为啥还要加一个wsum?
wsum(v, w::AbstractVector, [dim])
wsum!(R::AbstractArray, A::AbstractArray, w::AbstractVector, dim::Int; init::Bool=true)

# Statistics.mean 和 mean!
mean(A::AbstractArray, w::AbstractWeights[, dims::Int])
mean!(R::AbstractArray, A::AbstractArray, w::AbstractWeights[; dims=nothing])

进阶平均值计算

geomean(a)    # 几何平均
harmmean(a)   # 调和平均
genmean(a, p) # 幂平均, 当p = 0时, 就是几何平均
算术平均、几何平均、调和平均、平方平均和移动平均

  • 算术平均: 人人都知道的平均, 适用于数值型数据, 不太容易受随机因素影响(相比于中位数和众数), 但很容易受极值影响

  • 几何平均: 计算公式: \(G(x_1, x_2, ..., x_n) = \sqrt[n]{\prod_{i=1}^n X_i} \), 即连乘后开n次方, 适用于比率, 指数等的平均, 如增长率, 对数正态分布

  • 调和平均: 计算公式为: \( H = \cfrac{n}{\sum_{i=1}^n \cfrac{1}{x_i}} \), 即 倒数的算术平均再取倒数, 一般用于计算平均速率

以上三种平均统称毕达哥拉斯平均

矩/动差(moments)和积累量(cumulants)

  • 扩展了Statistics中的var,std以支持权重

  • mean_and_var,mean_and_std: 同时计算mean和var/std

  • skewness(v, [wv::AbstractWeights], m=mean(v)): 标准偏斜度

  • kurtosis(v, [wv::AbstractWeights], m=mean(v)): 峭度

  • moment(v, k, [wv::AbstractWeights], m=mean(v)): 中心矩

  • cumulant(v, k, [wv::AstractWeights], m=mean(v)): 积累量

中心矩, 积累量

中心矩\(μ_k\): 又叫主动差/中央动差(矩又叫做动差)

0阶中心矩 -> 1 1阶中心矩 -> 0 2阶 -> 方差 3阶 -> 偏度 4阶 -> 峰度

性质

  • 平移不变

  • n阶中心矩是n次齐次函数

  • 只有阶数小于等于3, 且两个相互独立的随机变量时, 其中心矩才具有加法性: \( \mu_n(X + Y) = \mu_n(X) + \mu_n(Y) \)

积累量:

:看维基吧

一阶积累量 -> 期望值 二阶积累量 -> 方差 三阶积累量 -> 三阶中心矩 四阶以上 -> 与中心矩不同

变异

  • span(x): 返回一个集合的范围: minimum(x):maximum(x)

  • variation(x, m=mean(x)): 变异系数 \(\frac{sd}{mean}\)

  • sem(x; mean=nothing): 标准误, 支持权重向量类型, FrequencyWeights和ProbablilityWeights的计算公式不同

  • mad(x; center=median(x), normalize=true)和mad!(语法一致, 结果替换x): 绝对中位差, 反映数据的波动情况。默认中心值为中位数, 默认normalize为true, 即在假设数据服从正态分布的前提下, 把原始mad值乘以1 / quantile(Normal(), 3/4) ≈ 1.4826从而估计标准差。

Z-score

zscore(X, [μ, σ])
# z-score(x_i) = (x_i - μ) / σ
# μ: mean(X); σ: std(X)

zscore!([Z], X, μ, σ)
# 如果提供了Z变量, 则将结果存到Z, 但是Z必须与X是同类型的
# 如果不提供Z, 则覆盖X

熵相关

entropy(p, [b])
# 计算一组概率值p的信息熵, 可以显式设置b用于scale: 1/log(b)

renyientropy(p, α)
# 计算向量p在顺序α下的瑞丽熵

crossentropy(p, q, [b])
# 计算p和q的交叉熵, 可以设置b用于scale

kldivergence(p, q, [b])
# 计算相对熵(K-L散度)
香农熵, 瑞丽熵, 交叉熵, 相对熵 TBC

分位数相关

percentile(x, p) # 返回x集合中p百分位数的值, 等同于: quantile(x, p / 100)
iqr(x) # 计算四分差, percentile(x, 75) - percentile(x, 25)
nquantile(x, n::Integer) # 返回能把x分成n等份的分位数集合
quantile(v, w::AbstractWeights, p) # 扩展Statistics的quantile函数, 支持权重
median(v, w::AbstractWeights) # 扩展median函数,支持权重

quantilerank(itr, value; method=:inc) #返回value在itr中的分位数排名,在[0,1]范围内的比值
# 设:
#    cl = itr中小于value的数目, count_less
#    ce = itr中等于value的数目, count_equal
#    n  = itr的长度
#    gs = itr中小于value的元素的最大值, greatest_smaller
#    sg = itr中大于value的元素的最小值, smallest_greater

# method可选参数和意义:
# :inc => 返回[0, 1]范围的值, 如果value ∈ itr, 返回cl/(n-1), 
#         否则按照Hyndman and Fan (1996)中定义7来计算, 与
#         Excel中的PERCENTRANK函数行为一致
# :exc => 返回(0, 1)范围的值, 如果value ∈ itr, 返回(cl+1)/(n+1),
#         否则按照Hyndman and Fan (1996)中定义6来计算, 
#         与Excel中的PERCENTRANK.EXC行为一致
# :compete => 如果value ∈ itr, 返回cl/(n-1), 否则返回(cl-1)/(n-1),
#         与MariaDB的PERCENT_RANK和dplyr的percent_rank行为一致
# :tied => 返回(cl + ce/2)/n, 与"mean"类型的SciPy的percentileofscore行为一致
# :strict => 返回 cl/n, 与"strict"类型的SciPy的percentileofscore行为一致
# :weak => 返回(cl + ce) / n, 与"weak"类型的SciPy的percentileofscore行为一致

percentilerank(itr, value; method=:inc) # 等同于 quantilerank(itr, value) * 100

众数

mode(a, [r]) # r => range
mode(a, w::AbstractWeights)
# 返回一组范围或者指定权重下的众数, 如果有多组结果, 只返回第一个

modes(a, [r])::Vector
modes(a, w)::Vector
# 返回所有众数结果向量

概括统计量

summarystats(a) # 统计a的长度, 缺失值, mean和[0,25,50,75,100]的百分位数:
julia> summarystats(1:100)
Summary Stats:
Length:         100
Missing Count:  0
Mean:           50.500000
Minimum:        1.000000
1st Quartile:   25.750000
Median:         50.500000
3rd Quartile:   75.250000
Maximum:        100.000000

describe(a) # 格式化输出summarystats(a)的结果 (我也没看出来好看在哪里啊)

信度分析

cronbachalpha(covmatrix::AbstractMatrix{<:Real})
# 根据协方差矩阵计算Cronbach's alpha(1951)
# 返回一个`ConbachAlpha`类型的对象, 包含如下特征:
# alpha: α系数
# dropped: 缺失值

# 一般认为 Cronbach's α > 0.7 表示条目之间的一致性较好

稳健性统计

trim(x::AbstractVector; prop=0.0, count=0)
trim!(x::AbstractVector; prop=0.0, count=0) # 原地修改x
# 返回一个迭代器, 包含除去上下prop或者count个数的元素的x集合
# mean(trim(x)) 计算修剪后的x的平均值
# trimvar(x; prop=0.0, count=0) 计算修剪后的x的var
trimvar(x; prop=0.0, count=0)

winsor(x::AbstractVector; prop=0.0, count=0)
winsor!(x::AbstractVector; prop=0.0, count=0)
# 温塞处理: 超出范围的数值变为范围内最后一个数值

偏差统计

conteq(a, b)
# 统计a和b中相等元素的个数, a和b长度要相同, 按照顺序对齐比对
# 类似于: [isequal(x...) for x in zip(a,b)] |> sum

countne(a, b)
# a和b中不相同的个数

sqL2dist(a, b) # 计算L2距离的方差
L2dist(a, b) # 计算L2距离
L1dist(a, b) # 计算L1距离
Linfdist(a, b) # 计算L∞距离, 又叫契比雪夫距离(Chebyshev distance)
gkldiv(a, b) # 计算KL散度
meanad(a, b) # 计算 平均绝对偏差, mean(abs, a - b)
maxad(a, b) # 计算最大绝对偏差, maxabs(a - b)
msd(a, b) # 计算均方差, mean(abs2, a - b)
rmsd(a, b; normalize=false) # 计算均方根差, sqrt(msd(a, b))
psnr(a, b, maxv) # 计算峰值信噪比

散布矩阵和协方差